library(tidyverse)
library(plotly)
library(leaflet)
library(ggmap)
single_color <- "#3182bd"
basic_theme <- theme(
  plot.title = element_text(size = 24, vjust = .5, hjust = .5),
  axis.title = element_text(size = 16, hjust = .5, vjust = .5),
  panel.background = element_rect(fill = "white"),
  legend.position = "none",
  panel.grid.major.y = element_line(color = "grey75", size = .1, linetype = "solid")
)
drivers <- feather::read_feather("drivers.feather")
incidents <- feather::read_feather("incidents.feather")
non_motorists <- feather::read_feather("non_motorists.feather")
combo <- feather::read_feather("combo.feather")
moco_north <- 39.35
moco_south <- 38.93
moco_east <- -76.89
moco_west <- -77.53
ns_center <- (moco_north + moco_south)/2
ew_center <- (moco_east + moco_west)/2
ns_size <- moco_north - moco_south
ew_size <- ns_size * 16/9
eee <- ew_center + (ew_size/2)
www <- ew_center - (ew_size/2)
pic <- get_stamenmap(bbox = c(left = www, bottom = moco_south, right = eee, top = moco_north), maptype = "toner-lines", zoom = 11) %>%
  ggmap(extent = "device",padding = 0) +
  theme(
    panel.background = element_blank(), panel.grid = element_blank(), 
    axis.title = element_blank(),axis.text = element_blank(),axis.ticks = element_blank()
    )
## Source : http://tile.stamen.com/toner-lines/11/582/780.png
## Source : http://tile.stamen.com/toner-lines/11/583/780.png
## Source : http://tile.stamen.com/toner-lines/11/584/780.png
## Source : http://tile.stamen.com/toner-lines/11/585/780.png
## Source : http://tile.stamen.com/toner-lines/11/586/780.png
## Source : http://tile.stamen.com/toner-lines/11/582/781.png
## Source : http://tile.stamen.com/toner-lines/11/583/781.png
## Source : http://tile.stamen.com/toner-lines/11/584/781.png
## Source : http://tile.stamen.com/toner-lines/11/585/781.png
## Source : http://tile.stamen.com/toner-lines/11/586/781.png
## Source : http://tile.stamen.com/toner-lines/11/582/782.png
## Source : http://tile.stamen.com/toner-lines/11/583/782.png
## Source : http://tile.stamen.com/toner-lines/11/584/782.png
## Source : http://tile.stamen.com/toner-lines/11/585/782.png
## Source : http://tile.stamen.com/toner-lines/11/586/782.png
## Source : http://tile.stamen.com/toner-lines/11/582/783.png
## Source : http://tile.stamen.com/toner-lines/11/583/783.png
## Source : http://tile.stamen.com/toner-lines/11/584/783.png
## Source : http://tile.stamen.com/toner-lines/11/585/783.png
## Source : http://tile.stamen.com/toner-lines/11/586/783.png
#ggsave("pres_background.png", plot = pic, width = 12.80, height = 7.20, units = "in")
moco_map <- get_stamenmap(bbox = c(left = www, bottom = moco_south, right = eee, top = moco_north), maptype = "terrain", zoom = 11)
## Source : http://tile.stamen.com/terrain/11/582/780.png
## Source : http://tile.stamen.com/terrain/11/583/780.png
## Source : http://tile.stamen.com/terrain/11/584/780.png
## Source : http://tile.stamen.com/terrain/11/585/780.png
## Source : http://tile.stamen.com/terrain/11/586/780.png
## Source : http://tile.stamen.com/terrain/11/582/781.png
## Source : http://tile.stamen.com/terrain/11/583/781.png
## Source : http://tile.stamen.com/terrain/11/584/781.png
## Source : http://tile.stamen.com/terrain/11/585/781.png
## Source : http://tile.stamen.com/terrain/11/586/781.png
## Source : http://tile.stamen.com/terrain/11/582/782.png
## Source : http://tile.stamen.com/terrain/11/583/782.png
## Source : http://tile.stamen.com/terrain/11/584/782.png
## Source : http://tile.stamen.com/terrain/11/585/782.png
## Source : http://tile.stamen.com/terrain/11/586/782.png
## Source : http://tile.stamen.com/terrain/11/582/783.png
## Source : http://tile.stamen.com/terrain/11/583/783.png
## Source : http://tile.stamen.com/terrain/11/584/783.png
## Source : http://tile.stamen.com/terrain/11/585/783.png
## Source : http://tile.stamen.com/terrain/11/586/783.png
temp_data <- incidents %>%
  select(acrs_report_type_fatal_crash, latitude, longitude)
ggmap(moco_map) +
  geom_point(data = temp_data %>% filter(acrs_report_type_fatal_crash == 0), aes(x = longitude, y = latitude), size = .2, color = "black") +
  geom_point(data = temp_data %>% filter(acrs_report_type_fatal_crash == 1), aes(x = longitude, y = latitude), size = 1.5, color = "red") +
  theme(
    legend.position = "none",
    axis.text = element_blank(), axis.ticks = element_blank(), axis.title = element_blank(),
    panel.background = element_blank(), panel.grid = element_blank(),
    plot.title = element_text(size = 20 , hjust = .5)
  )
## Warning: Removed 105 rows containing missing values (geom_point).

temp <- incidents %>%
  filter(!not_in_county) %>%
  mutate(non_motorist_involved = !(related_non_motorist == "N/A"))
temp1 <- temp %>% filter(non_motorist_involved)
temp2 <- temp %>% filter(!non_motorist_involved)
leaflet() %>%
  addTiles() %>%
  addCircleMarkers(
    data = temp2,
    lng = ~longitude, lat = ~latitude,
    radius = 3, stroke = TRUE, weight = 2, opacity = 1, color = "blue", group = "Vehicles"
  ) %>%
    addCircleMarkers(
    data = temp1,
    lng = ~longitude, lat = ~latitude,
    radius = 3, stroke = TRUE, weight = 2, opacity = 1, color = "red", group = "Non Motorists"
  ) %>%
  addLayersControl(overlayGroups = c("Non Motorists", "Vehicles"))
pal <- colorFactor(c("red","orange", "black"), domain = c("Driver", "Non Motorist", "Other"))
temp <- combo  %>%
  filter(acrs_report_type_fatal_crash == 1) %>%
  mutate(
    fatal_type = case_when(
      nm_injury_severity_fatal_injury > 0 ~ "Non Motorist",
      d_injury_severity_fatal_injury > 0 ~ "Driver",
      TRUE ~ "Other"
    )
    )
leaflet() %>%
  addTiles() %>%
  addCircleMarkers(
    data = temp,
    lng = ~longitude, lat = ~latitude,
    radius = 3, stroke = TRUE, weight = 2, opacity = 1, color = ~pal(fatal_type), group = ~fatal_type,
    popup = ~paste(crash_date_time, road_name, sep = "<br/>")
  ) %>%
  addLayersControl(overlayGroups = c("Driver","Non Motorist","Other"))
ordered_injury <- drivers %>%
  count(injury_severity) %>%
  arrange(desc(n)) %>%
  pull(injury_severity)
speed_limit_plot <- drivers %>%
  mutate(
    speed_limit = factor(c("LOW", "MEDIUM", "HIGH")[1 + findInterval(speed_limit, c(31, 46))], levels = c("LOW", "MEDIUM", "HIGH")),
    injury_severity = factor(injury_severity, ordered_injury)
  ) %>%
  group_by(speed_limit, injury_severity) %>%
  tally() %>%
  rename(`Speed Limit` = speed_limit, `Injury Severity` = injury_severity) %>%
  ggplot() +
  geom_bar(aes(x = `Speed Limit`, fill = `Injury Severity`,  weight = n), position = "fill") +
  basic_theme +
  theme(legend.position = "top") +
  ggtitle("Injuries and Speed Limits") +
  ylab("Proportion")
ggplotly(speed_limit_plot)
ordered_injury <- drivers %>%
  count(injury_severity) %>%
  arrange(desc(n)) %>%
  pull(injury_severity)
speed_limit_plot2 <- drivers %>%
  mutate(
    speed_limit = factor(as.character(speed_limit), levels = as.character(seq.int(from = 0, to = 70, by = 5))),
    injury_severity = factor(injury_severity, ordered_injury)) %>%
  group_by(speed_limit, injury_severity) %>%
  tally() %>%
  rename(`Speed Limit` = speed_limit, `Injury Severity` = injury_severity) %>%
  ggplot() +
  geom_bar(aes(x = `Speed Limit`, fill = `Injury Severity`,  weight = n), position = "fill") +
  basic_theme +
  theme(legend.position = "top") +
  ggtitle("Injuries and Speed Limits") +
  ylab("Proportion")
ggplotly(speed_limit_plot2) %>% layout(legend = list(orientation = 'h'))